16LO

Technical details

Show code
library(GeoPressureR)
library(leaflet)
library(leaflet.extras)
library(raster)
library(dplyr)
library(ggplot2)
library(kableExtra)
library(plotly)
library(GeoLocTools)
setupGeolocation()
knitr::opts_chunk$set(echo = FALSE)
load(paste0("../data/1_pressure/", params$gdl_id, "_pressure_prob.Rdata"))
load(paste0("../data/2_light/", params$gdl_id, "_light_prob.Rdata"))
load(paste0("../data/3_static/", params$gdl_id, "_static_prob.Rdata"))
load(paste0("../data/4_basic_graph/", params$gdl_id, "_basic_graph.Rdata"))

Settings used

All the results produced here are generated with (1) the raw geolocator data, (2) the labeled files of pressure and light and (3) the parameters listed below.

Show code
kable(gpr) %>% scroll_box(width = "100%")
gdl_id crop_start crop_end thr_dur extent_N extent_W extent_S extent_E map_scale map_max_sample map_margin prob_map_s prob_map_thr shift_k kernel_adjust calib_lon calib_lat calib_1_start calib_1_end calib_2_start calib_2_end calib_2_lon calib_2_lat prob_light_w thr_prob_percentile thr_gs RingNo scientific_name common_name mass wing_span Color
16LO 2017-01-10 2017-12-23 6 17 9 -25 38 5 300 30 1.2 0.9 0 1.4 28.77132 -22.72004 2017-01-10 2017-03-26 2017-11-15 2017-12-23 NA NA 0.1 0.95 120 NA Halcyon senegaloides Woodland Kingfisher NA NA #FF70A6

Pressure timeserie

The labeling of pressure data is illustrated with this figure. The black dots indicates the pressure datapoint not considered in the matching. Each stationary period is illustrated by a different colored line.

Show code
pressure_na <- pam$pressure %>%
  mutate(obs = ifelse(isoutliar | sta_id == 0, NA, obs))
p <- ggplot() +
  geom_line(data = pam$pressure, aes(x = date, y = obs), colour = "grey") +
  geom_point(data = subset(pam$pressure, isoutliar), aes(x = date, y = obs), colour = "black") +
  # geom_line(data = pressure_na, aes(x = date, y = obs, color = factor(sta_id)), size = 0.5) +
  geom_line(data = do.call("rbind", shortest_path_timeserie) %>% filter(sta_id > 0), aes(x = date, y = pressure0, col = factor(sta_id))) +
  theme_bw() +
  scale_colour_manual(values = pam$sta$col) +
  scale_y_continuous(name = "Pressure(hPa)")

ggplotly(p, dynamicTicks = T) %>% layout(showlegend = F)

Pressure calibration

Show code
sp_pressure = do.call("rbind", shortest_path_timeserie) %>% filter(sta_id > 0)

sta_plot <- which(difftime(pam$sta$end,pam$sta$start,unit="days")>3)

par(mfrow=c(2,3))
for (i in seq_len(length(sta_plot))){
  i_s = sta_plot[i]
  pressure_s = pam$pressure %>% 
    filter(sta_id==i_s & !isoutliar)
  
    err <- pressure_s %>% left_join(sp_pressure, by="date") %>% 
      mutate(
        err = obs-pressure-mean(obs-pressure)
      ) %>% .$err
    
    hist(err, freq = F, main = paste0("sta_id=",i_s, " | ",nrow(pressure_s)," dtpts | std=",round(sd(err),2)))
   xfit <- seq(min(err), max(err), length = 40) 
    yfit <- dnorm(xfit, mean = mean(err), sd = sd(err)) 
    lines(xfit, yfit, col = "red")
}

Light

Show code
raw_geolight <- pam$light %>%
  transmute(
    Date = date,
    Light = obs
  )
lightImage(tagdata = raw_geolight, offset = 0)
tsimagePoints(twl$twilight,
  offset = 0, pch = 16, cex = 1.2,
  col = ifelse(twl$deleted, "grey20", ifelse(twl$rise, "firebrick", "cornflowerblue"))
)
abline(v = gpr$calib_2_start, lty = 1, col = "firebrick", lwd = 1.5)
abline(v = gpr$calib_1_start, lty = 1, col = "firebrick", lwd = 1.5)
abline(v = gpr$calib_2_end, lty = 2, col = "firebrick", lwd = 1.5)
abline(v = gpr$calib_1_end, lty = 2, col = "firebrick", lwd = 1.5)

Show code
hist(z, freq = F)
lines(fit_z, col = "red")

The probability map resulting from light data alone can be seen below.

Show code
li_s <- list()
l <- leaflet(width = "100%") %>%
  addProviderTiles(providers$Stamen.TerrainBackground) %>%
  addFullscreenControl()
for (i_r in seq_len(length(light_prob))) {
  i_s <- metadata(light_prob[[i_r]])$sta_id
  info <- pam$sta[pam$sta$sta_id == i_s, ]
  info_str <- paste0(i_s, " | ", info$start, "->", info$end)
  li_s <- append(li_s, info_str)
  l <- l %>% addRasterImage(light_prob[[i_r]], opacity = 0.8, colors = "OrRd", group = info_str)
}
l %>%
  addCircles(lng = gpr$calib_lon, lat = gpr$calib_lat, color = "black", opacity = 1) %>%
  addLayersControl(
    overlayGroups = li_s,
    options = layersControlOptions(collapsed = FALSE)
  ) %>%
  hideGroup(tail(li_s, length(li_s) - 1))

Light vs Pressure

We can compare light and pressure location at long stationary stopover (>5 days). By assuming the best match of the pressure to be the truth, we can plot the histogram of the zenith angle and compare to the fit of kernel density at the calibration site.

Show code
 raw_geolight <- pam$light %>%
    transmute(
      Date = date,
      Light = obs
    )
 dur <- unlist(lapply(pressure_prob, function(x) difftime(metadata(x)$temporal_extent[2],metadata(x)$temporal_extent[1], units = "days" )))
  long_id <- which(dur>5)

par(mfrow = c(2, 3))
for (i_s in long_id){
  twl_fl <- twl %>%
    filter(!deleted) %>%
    filter(twilight>shortest_path_timeserie[[i_s]]$date[1] & twilight<tail(shortest_path_timeserie[[i_s]]$date,1))
  sun <-  solar(twl_fl$twilight)
  z_i <- refracted(zenith(sun, shortest_path_timeserie[[i_s]]$lon[1], shortest_path_timeserie[[i_s]]$lat[1]))
  hist(z_i, freq = F, main = paste0("sta_id=",i_s, " | ",nrow(twl_fl),"twls"))
  lines(fit_z, col = "red")
  xlab("Zenith angle")
}

Similarly, we can plot the line of sunrise/sunset at the best match of pressure (yellow line) and compare to the raw and labeled light data.

Show code
  lightImage(
    tagdata = raw_geolight,
    offset = gpr$shift_k / 60 / 60
  )
  tsimagePoints(twl$twilight,
                offset = gpr$shift_k / 60 / 60, pch = 16, cex = 1.2,
                col = ifelse(twl$deleted, "grey20", ifelse(twl$rise, "firebrick", "cornflowerblue"))
  )
  for (ts in shortest_path_timeserie){
    if (!is.null(ts)){
      twl_fl <- twl %>%
      filter(twilight>ts$date[1] & twilight<tail(ts$date,1))
      if (nrow(twl_fl)>0){
      tsimageDeploymentLines(twl_fl$twilight,
                             lon = ts$lon[1], ts$lat[1],
                             offset = gpr$shift_k / 60 / 60, lwd = 3,col = adjustcolor("orange", alpha.f = 0.5))
        
      }
    }
  }

GeoPressureViz

To visualize the path on GeoPressureViz, you will need to also load the pressure and light probability map and align them first with the code below.

Show code
sta_marginal <- unlist(lapply(static_prob_marginal, function(x) raster::metadata(x)$sta_id))
sta_pres <- unlist(lapply(pressure_prob, function(x) raster::metadata(x)$sta_id))
sta_light <- unlist(lapply(light_prob, function(x) raster::metadata(x)$sta_id))
pressure_prob <- pressure_prob[sta_pres %in% sta_marginal]
light_prob <- light_prob[sta_light %in% sta_marginal]

The code below will open with the shortest path computed with the graph approach.

Show code
geopressureviz <- list(
  pam_data = pam,
  static_prob = static_prob,
  static_prob_marginal = static_prob_marginal,
  pressure_prob = pressure_prob,
  light_prob = light_prob,
  pressure_timeserie = shortest_path_timeserie
)
save(geopressureviz, file = "~/geopressureviz.RData")

shiny::runApp(system.file("geopressureviz", package = "GeoPressureR"),
  launch.browser = getOption("browser")
)

Stationay period information

Show code
pam$sta %>% mutate(duration = difftime(end,start,units="days")) %>% kable()
start end sta_id col duration
2017-01-10 00:00:00 2017-03-25 17:20:00 1 #1B9E77 74.7222222 days
2017-03-26 03:45:00 2017-03-26 18:30:00 2 #D95F02 0.6145833 days
2017-03-26 20:25:00 2017-03-28 20:30:00 3 #7570B3 2.0034722 days
2017-03-29 03:20:00 2017-03-29 19:55:00 4 #E7298A 0.6909722 days
2017-03-29 20:45:00 2017-03-30 21:45:00 5 #66A61E 1.0416667 days
2017-03-31 00:25:00 2017-03-31 21:40:00 6 #E6AB02 0.8854167 days
2017-04-01 00:00:00 2017-04-01 17:50:00 7 #A6761D 0.7430556 days
2017-04-01 20:35:00 2017-04-02 00:20:00 8 #666666 0.1562500 days
2017-04-02 01:00:00 2017-04-02 19:50:00 9 #1B9E77 0.7847222 days
2017-04-02 21:45:00 2017-04-04 01:45:00 10 #D95F02 1.1666667 days
2017-04-04 03:00:00 2017-04-05 01:15:00 11 #7570B3 0.9270833 days
2017-04-05 02:15:00 2017-04-06 23:40:00 12 #E7298A 1.8923611 days
2017-04-07 01:25:00 2017-04-08 23:00:00 13 #66A61E 1.8993056 days
2017-04-09 00:15:00 2017-04-09 22:15:00 14 #E6AB02 0.9166667 days
2017-04-09 22:55:00 2017-04-15 19:45:00 15 #A6761D 5.8680556 days
2017-04-16 01:20:00 2017-04-21 17:55:00 16 #666666 5.6909722 days
2017-04-22 03:45:00 2017-04-22 21:35:00 17 #1B9E77 0.7430556 days
2017-04-23 03:45:00 2017-04-23 18:30:00 18 #D95F02 0.6145833 days
2017-04-24 02:05:00 2017-04-24 20:40:00 19 #7570B3 0.7743056 days
2017-04-24 21:20:00 2017-04-25 21:55:00 20 #E7298A 1.0243056 days
2017-04-26 01:15:00 2017-04-27 00:40:00 21 #66A61E 0.9756944 days
2017-04-27 03:25:00 2017-04-27 18:55:00 22 #E6AB02 0.6458333 days
2017-04-27 22:25:00 2017-04-29 00:50:00 23 #A6761D 1.1006944 days
2017-04-29 02:40:00 2017-04-29 22:30:00 24 #666666 0.8263889 days
2017-04-30 00:30:00 2017-05-30 20:20:00 25 #1B9E77 30.8263889 days
2017-05-30 22:35:00 2017-06-01 20:25:00 26 #D95F02 1.9097222 days
2017-06-01 22:20:00 2017-10-16 17:10:00 27 #7570B3 136.7847222 days
2017-10-17 03:25:00 2017-10-17 16:30:00 28 #E7298A 0.5451389 days
2017-10-18 03:25:00 2017-10-18 16:25:00 29 #66A61E 0.5416667 days
2017-10-19 01:40:00 2017-10-19 18:40:00 30 #E6AB02 0.7083333 days
2017-10-19 22:40:00 2017-10-20 03:05:00 31 #A6761D 0.1840278 days
2017-10-20 03:25:00 2017-10-20 18:30:00 32 #666666 0.6284722 days
2017-10-21 03:25:00 2017-10-21 19:15:00 33 #1B9E77 0.6597222 days
2017-10-21 23:40:00 2017-10-22 21:35:00 34 #D95F02 0.9131944 days
2017-10-23 00:30:00 2017-10-24 00:05:00 35 #7570B3 0.9826389 days
2017-10-24 02:25:00 2017-10-25 00:20:00 36 #E7298A 0.9131944 days
2017-10-25 02:25:00 2017-10-26 18:40:00 37 #66A61E 1.6770833 days
2017-10-26 21:00:00 2017-10-27 21:25:00 38 #E6AB02 1.0173611 days
2017-10-28 01:20:00 2017-10-31 18:45:00 39 #A6761D 3.7256944 days
2017-10-31 23:30:00 2017-11-01 17:30:00 40 #666666 0.7500000 days
2017-11-01 23:55:00 2017-11-02 19:45:00 41 #1B9E77 0.8263889 days
2017-11-02 22:25:00 2017-11-02 23:55:00 42 #D95F02 0.0625000 days
2017-11-03 02:10:00 2017-11-04 01:40:00 43 #7570B3 0.9791667 days
2017-11-04 02:50:00 2017-11-07 23:10:00 44 #E7298A 3.8472222 days
2017-11-08 01:00:00 2017-11-08 23:10:00 45 #66A61E 0.9236111 days
2017-11-09 00:15:00 2017-11-09 19:40:00 46 #E6AB02 0.8090278 days
2017-11-09 21:00:00 2017-11-10 21:25:00 47 #A6761D 1.0173611 days
2017-11-11 00:10:00 2017-11-12 17:50:00 48 #666666 1.7361111 days
2017-11-12 22:15:00 2017-11-13 20:20:00 49 #1B9E77 0.9201389 days
2017-11-14 01:30:00 2017-11-14 19:00:00 50 #D95F02 0.7291667 days
2017-11-15 00:15:00 2017-12-22 23:55:00 51 #7570B3 37.9861111 days